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Using modern methods of reactor physics we have performed full-scale calculations of the natural 
reactor Oklo. For reliability we have used recent version of two Monte Carlo codes: Russian code 
MCU REA and world wide known code MCNP (USA) . Both codes produce similar results. We have 
constructed a computer model of the reactor Oklo zone RZ2 which takes into account all details 
of design and composition. The calculations were performed for three fresh cores with different 
uranium contents. Multiplication factors, reactivities and neutron fluxes were calculated. We have 
estimated also the temperature and void effects for the fresh core. As would be expected, we have 
1^-^ ' found for the fresh core a significant difference between reactor and Maxwell spectra, which was 

I used before for averaging cross sections in the Oklo reactor. The averaged cross section of Sm 

. and its dependence on the shift of resonance position (due to variation of fundamental constants) 

' are significantly different from previous results. 

Contrary to results of some previous papers we find no evidence for the change of the fine structure 
constant in the past and obtain new, most accurate limits on its variation with time: 
—4 • 10~^^year ~^ < a/a < 3 • 10~^^ year~^. A further improvement in the accuracy of the limits 
can be achieved by taking account of the core burnup. These calculations are in progress. 
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OO ■ I. INTRODUCTION 

' The discovery of the natural nuclear reactor in Gabon (West Africa) was possibly one of the most momentous events 
in reactor physics since in 1942 Enrico Fermi with his team achieved an artificial self-sustained fission chain reaction. 
CD ' Soon after the discovery of the ancient natural reactor Oklo in Gabon (West Africa) Q, |^ |1] , one the authors of 
, the present paper (Yu.P.) and his postgraduate student A.I. Shlyakhter reahzed that the "Oklo phenomenon" could 
O • be used to find the most precise limits on possible changes of fundamental constants. At that time they considered 
I ] probabilistic predictions of unknown absorption cross sections based on static nuclear properties. Near the neutron 
binding energy (_B„ = 6 — 8 MeV) the resonances form a fence with a mean separation of tens of electron volts. The 
(-H I magnitude of the cross section depends on the proximity of the energy Eth of the thermal neutron to the nearest 
• • . resonance. If the energy Eth = 25 meV falls directly on a resonance, then the cross section increases to as much 
10^ — 10^ b If the cross sections have changed with time, then the entire fence of resonances as a whole has 

. shifted by a small amount AEr. This shift can be established most accurately from the change of the cross section 
' of strong absorbers (for instance of ^62 ^^s)- The estimate of the shift AEr on the basis of experimental data for the 
5^ , Oklo reactor {AEr < 3 • 10^^'' eV/year) 0,|1,ISE1 allowed us to get the most accurate estimate of a possible limit of 
the rate of change of fundamental constants. This estimate has remained the most accurate one for 20 years. In 1996 
a paper by Damour and Dyson was published in which the authors checked and confirmed the results of Ref . |3, ISj, Q . 
Damour and Dyson were the first to calculate the dependence of the capture cross section on the temperature Tc 
of the core: ct^ sm (Tc) [Hill. In 2000 Fujn et al. published a paper in which the authors significantly reduced 
the experimental error of the cross section (Jr,Sm{Tc) uM- both papers the authors averaged the samarium cross 
section with a Maxwell velocity spectrum over a wide interval of the core temperature Tc- 

After the publication of Refs. [l3)ll3li ^^^^ the authors of the present paper (Yu. P.) realized that the limit on 
the change of the cross section can be significantly improved at least in two directions: 

1. Instead of a Maxwell distribution, the samarium cross section should be averaged with the spectrum of the Oklo 
reactor which contains the tail of Fermi spectrum of slowing down cpithcrmal neutrons. 

2. The range of admissible core temperatures Tc can be significantly reduced, assuming that Tc is the equilibrium 
temperature at which the effective multiplication factor K^ff of the reactor is equal to one. On account of the 
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negative power coefficient (void + temperature) such a state of tfie reactor will be maintained for a long time 
until the burn- up results in a reduction of the reactivity excess and hence of Tq- Since ^gf^™ burns up about 
100 times faster than ^gfU, the core will contain only that amount of samarium that was generated immediately 
before the reactor shut down. Therefore one needs to know the reactivity excess and Tc at the end of the cycle. 

To solve this problem one must use modern neutron-physical and thermo-hydrodynamical methods of reactor 
calculations. We have built up a complete computer model of the Oklo reactor core RZ2 and established its material 
composition. We have chosen three variants of its initial composition in order to estimate its effect on the spread 
of results. To increase the reliability of the results wc have used modern versions of two Monte Carlo codes. One 
of them, which has been developed at the Kurchatov Institute, is the licensed Russian code MCU-REA with the 
library DLC/MCUDAT-2.2 of nuclear data Q; the other one is the well known international code MCNP4C with 
library ENDF/B-VI (13 • Both codes give similar results. We have calculated the multiplication factors, reactivity and 
neutron flux for the fresh cores, and the void and temperature effects As expected, the reactor spectrum differs 
strongly from a Maxwell distribution (see below). The cross section ar,sm{Tci AE'r), averaged with this distribution, 
is significantly different from the cross section averaged with a Maxwell distribution We use our result for the 
averaged cross section to estimate the position of resonances at the time of Oklo reactor activity. This allows us to 
obtain the most accurate limits on the change of the fine structure constant in the past. 

The paper is organized as follows. In section I we describe briefly the history of the discovery of the natural Oklo 
reactor and itemize the main parameters of its cores. We consider mainly the core RZ2 . We describe in detail the 
neutronics of this core calculated by modern Monte Carlo codes. However simple semianalytical considerations are 
also useful to clarify the picture. We consider the power effect which is a sum of the temperature and void effects. At 
the end of the section we discuss the computational difficulties in the calculations of the unusually large core RZ2 
and demonstrate that Monte Carlo methods are, in general, inadequate for the calculations of core burn-up. 

The main result of Section I is the neutron spectrum in the fresh core. In Section II we apply this spectrum to 
obtain the averaged cross section of ^|2Sm in the past. We begin this Section with an explanation of the way of 
obtaining precise limits on the variation of fundamental constants using the available Oklo reactor data. We describe 
different approaches to the problem and relate the variation of the constants to the change in the averaged cross 
sections for thermal neutrons. Using our value for the cross section of ^62^™ obtain limits on the variation of the 
fine structure constant which is the best available at the moment. At the end of this Section we compare our result 
with the results obtained in other papers and discuss possible reasons of differences. 

II. NEUTRONICS OF THE FRESH CORE 
A. History of the discovery and parameters of the Oklo reactor 

1. History of the discovery of the natural reactor 

The first physicist to say in May or June of 1941 that a nuclear chain reaction could have been more easily realized 
a billion years ago was Yakov Borisovich Zcldovich At that time he was considering the possibility of getting 
a fission chain reaction in a homogeneous mixture of natural uranium with ordinary water. His calculations (with 
Yu.B. Khariton) showed that this could be achieved with an approximately two- fold enrichment of natural uranium 
|l9l l20|. A billion years ago the relative concentration of the light uranium isotope was significantly higher, and a 
chain reaction was possible in a mixture of natural uranium and water. "Yakov Borisovich said nothing about the 
possibility of a natural reactor, but his thoughts directly lead us to the natural reactor discovered in Gabon in 1972" 
reminisced I.I. Gurevich Later, in 1957, G. Whetherill and M. Inghram arrived at the same conclusion [2ll[2^ . 
Going from the present concentration of uranium in pitchblende, they concluded that about two billion years ago, 
when the proportion of ^gfU exceeded 3%, conditions could be close to critical. Three years later, P. Kuroda [23Ll2^ 
showed that, if in the distant past there was water present in such deposits, then the neutron multiplication factor 
(Koo) for an infinite medium could exceed unity and a spontaneous chain reaction could arise. But before 1972 no 
trace of a natural reactor has been found. On the 7th of June 1972, during a routine mass-spectroscopic analysis 
in the French Pierrelatte factory that produced enriched fuel, H. Bouzigues 0, Q noticed that the initial uranium 
hexafluoride contains Cs = 0.717% of ^g2U atoms instead of the 0.720%, which is the usual concentration in terrestrial 
rock, meteorites and lunar samples. The French Atomic Energy Authority (CEA) began an investigation into this 
anomaly. The phenomenon was named " Oklo phenomenon" . The results of this research were published in the 
proceedings of two IAEA symposia 0, . The simplest hypothesis of a contamination of the uranium by depleted 
tails of the separation process was checked and shown to be wrong. Over a large number of steps of the production 
process, the anomaly was traced to the Munana factory near Franceville (Gabon) where the ore was enriched. The ore 
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with a mean uranium concentration of (0.4 — 0.5)% got delivered there from the Oklo deposit. The isotope analysis 
of the uranium- rich samples showed a significant depletion of the ^||U isotope and also a departure from the natural 
distribution of those rare earth isotopes, which are known as fission products 0,0,11^ 1131 ■ This served as a proof of 
the existence in the distant past of a spontaneous chain reaction. It had taken less than three months to produce this 
proof. A retrospective analysis of documents and samples of the Munanc enrichment factory showed that in 1970-72 
ore was delivered for processing that contained at times up to 20% of uranium depleted to 0.64% of isotope ^gfU [28l |. 
Considering that the ore was mixed during mining, the uranium concentration could be even higher in some samples, 
and the depletion even stronger. Altogether more than 700 tons of depleted uranium has been mined that had taken 
part in the chain reaction. The deficit of ^gfU (that had not been noticed at first) was about 200 kg. By agreement 
with the Government of Gabon, the uranium ore production company of Franceville (COMUF) agreed to halt mining 
in the region of the natural reactor. A Franco-Gabon group headed by R. Naudet began a systematic study of the 
Oklo phenomenon. Numerous samples, obtained by boring, were sent for analysis to various laboratories around the 
world. They allowed a reconstruction of the functioning of the reactor in the Precambrian epoch. 



2. Geological history of the Oklo deposit 

As was shown by the U/Pb analysis, the Oklo deposit with a uranium concentration of about 0.5% in the sediment 
layer was formed about 2 • 10^ years ago [2^|23j0|- During this epoch an important biological process was taking 
place: the transition from prokaryotes, i.e. cells without nucleus, to more complex unicellular forms containing a 
nucleus - eucaryotes. The eucaryotes began to absorb carbon oxide and hence saturate the atmosphere with oxygen. 
Under the influence of oxygen, the uranium oxides began transforming into forms containing more oxygen, which are 
soluble in water. Rains have washed them into an ancient river, forming in its mouth a sandstone sediment, rich in 
uranium, of 4 to 10 meters thickness and a width of 600 to 900 meters [sj. The heavier uranium particles settled 
more quickly to the ground of the nearly stagnant water of the river delta. As a result the sandstone layer got enriched 
with uranium up to 0.5% (as in an enrichment factory). After its formation, the uranium-rich layer, that was resting 
on a basalt bed, was covered by sediments and sank to a depth of 4 kilometers. The pressure on this layer was 100 
MPals^- Under this pressure the layer got fractured and ground water entered the clefts. Under the action of the 
filtered water that was subjected to a high pressure, and as a result of not completely understood processes, lenses 
formed with a very high uranium concentration (up to 20 — 60% in the ore) with a width of 10 to 20 meters and of the 
order of 1 meter thickness [s^l. The chain reaction took place in these lenses. After the end of the chain reaction the 
deposit was raised to the surface by complicated tectonic processes and became accessible for mining. Within tens of 
meters six centers of reactions were found immediately, and altogether the remains of 17 cores were found 34] . 

The age Tq of the reactor was determined from the total number of ^||U nuclei burnt up in the past, Nr,i,(^d) , and 
the number of nuclei existing today, A^5(To)(here is the density of ^gfU and d is the duration of the chain reaction). 
For such a way of determining Tq it is necessary to know the number of ^gfPu nuclei formed as a result of neutron 
capture by ^gfC/ and decayed to ^gfU, and the fluence 5* = ($ being the neutron flux). Another independent 
method consists of the determination of the amount of lead formed as a result of the decay of ^gf U, assuming that 
it did not occur in such a quantity in the initial deposit Both methods yield Tq — 1.81(5) • 10^ years [Tol Issj . 
Below we assume in our calculations the value of Tq = 1.8 • 10^ years. 

The duration of the work of the reactor can be established from the amount of ^g|Pu formed. One can separate 
the decayed ^gfPu from the decayed ^gfU using the different relative yields of Nd isotopes: 

= i50Nd/(i''3Nd-H "^Nd) = 0.1175 for ^IfPu and = 0.0566 for ^sfu ^j. However this comparison is masked 
by the fission of ^gfU by fast neutrons: (5^^ — 0.1336. Taking account of this contribution one arrives at an estimate 
of d '--^ 0.6 million years [33 . This was the value we adopted in our calculations. 

The total energy yield of the reactor has been estimated to be 1.5 • 10^ MWa [s^. Such a fission energy is obtained 
by two blocks of the Leningrad atomic power station with a hundred percent load in 2.3 years. Assuming a mean 
duration of d 6 • 10^ a for the work of the reactor one gets a mean power output of only Pp — 25 kW. 



B. Composition and size of the Oklo RZ2 reactor 

The cores of the Oklo reactor have been numbered. The most complete data are available for core RZ2 . This core 
of the Oklo reactor is of the shape of an irregular rectangular plate that lies on a basalt bed at an angle of 45°. The 
thickness of the plate is = 1 m, its width is 6 = 11 — 12 m, and its length is Z = 19 — 20 m (see Fig. 8a in Refs.[3^ 
and [l^). Thus the volume of the RZ2 core is about 240 m^. Since in the case of large longitudinal and transverse 
sizes the shape of the reactor is not essential, we have assumed as a reactor model a flat cylinder of height H = 1 m 
and radius R which is determined by the core burn-up. The energy yield is Ppd — 1.5 • 10^ MWa ~ 5.48 • 10^ MWd. 
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At a consumption of ^gfU of g = 1.3 g/MWd|40l|. the total amount of burnt up fissile matter is 

AMfc = gPpd = 7.12 tons. 



(1) 



Taking into account that half of the burnt up g^U isotope is replenished from the decay of the produced 94PU, we 
find the original mass of the burnt up ^gf U: 

AM5 = 4.75 tons. 

In case of a uniform burn-up, the average density of the burnt up '^g2U is 

A75(rf) 



AM. 



where R is given in meters. The relative average initial burn-up is 



1.51 3 

g/cm , 



(2) 



(3) 



A75(d) 1.51 g/cm^ 



75,^(0) 75,^(0)i?f 



and Ri 



1.51 g/cm^ 



1/2 



(4) 



75,. (0)^5 (rf) 

Processing the data of Table 2 from Ref.jS^ gives a value of ^^(To) ~ 50% for the present-day average over the core. 
In the past it was 1.355 times smaller (see below) on account of the higher concentration 75,^(0) of uranium, i.e. 
^5(0) = 36.9%. Thus the radius is given by the following formula: 



Ri 



4.09 g/cm^ 
75,,: (0) 



1/2 



(5) 



The approximate composition of the rock in core RZ2 is shown in Table 1 of Ref.|0. On the basis of these data one 
can calculate the elemental composition of the ore by weight (see the penultimate column of Tabl^. For comparison 
we show in the last column the composition by weight from the book of Yu. A. Shukolyukov (Table 2.1 of Ref. Elj). 
which was based on early data of R. Naudet ,4^]. These values coincide within 10%. Since all those elements are 
relatively weakly absorbing, such differences practically do not play any role. 



TABLE I: Present-day composition of the empty rock [l^ 
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SiOa 


43.00 





15.999 


44.04 


44 
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AI2O3 


25.73 


Si 


28.086 


20.10 


20 
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FeO 


14.53 


Al 


26.982 


13.62 


16 
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Fe203 


4.47 


Fe 


55.847 


14.42 


11 


5 


MgO 


10.43 


Mg 


24.305 


6.30 


4 


6 


K2O 


1.84 


K 


39.098 


1.53 


2 


7 


Sum 


100 






100 


97 



It is much more important to know the amount of uranium and of water at the beginning of the work of the 
reactor. The connection between the uranium content in the core (in %) and the density of the dehydrated core has 
been measured experimentally in Ref. [43l | (Fig^ . The content of uranium in the core varies greatly between different 
samples. To determine the influence of the uranium content on the reactor parameters we have chosen three initial 
values for the density of uranium in the dehydrated ore: lc/,i(T'o) =35, 45 and 55%, taken to be constant over the 
reactor. The value of the ore density ji{To) that corresponds to yt/,i(To) is shown in the second row of Table ITU fsee 
Figffll- In the fifth row of Table ITU we show the density of the empty rock. The density of water in the reactor is 
0.3 — 0.5 g/cni"^ 34] . This water consists of bound (crystalline) and unbound water which evaporates after 100°C. In 
our reactor model we assumed a total density of water of 7/^20 — 0.355 g/cm^, of which 0.155 g/cm^ was taken for 
the density of unbound water. Assuming a porosity of about 6% one can take for the density of dry ore with water 
the same value as for the dehydrate ore. 

The density of 2g|U of the fresh core RZ2 at the epoch of the formation of the reactor was 



78(0) = 78(ro)(l - C5)exp(-HTo/r8), 



(6) 
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FIG. 1: Present-day abundance of uranium (in volume percent of dry ore) in dependence on the dry ore density 7u in the 
Oklo reactor |4^. Three initial variants of the Yu,i{To) = 35%, 45%, 55%U in dry ore are shown. 



where the lifetime of ^gfU is rg = 6.45 • 10 year. The value of 7^(0) increases on account of the decay of uranium 
into lead. The density of ^gfU (7(7,5(0), g/cm'^) in the fresh core (T5 — 1.015 • 10^ year) is 

7t/,5(0) = 7c/,5(ro)C5exp(+ro/r5). (7) 

The values of 75,1 (0) and 78,i(0) are shown in Table UTI Also in the table are the calculated values of the densities of 
uranium 7c/,i(0) and of the empty rock (without Pb), and the new fraction of uranium in the dry ore at the beginning 



TABLE II: Present-day and initial composition of the ore in core RZ2 
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Present-day fraction of U in the 
dry ore, % 


YuATo) 


35 


45 


55 


2 


Present-day density of the dry 
ore, g/cm^ (Fig0 


7i(ro) 


2.84 


3.29 


3.82 
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Present-day density of U in the 
dry ore, g/cm'' 


■yuATo) 


0.994 


1.481 


2.101 
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Present-day density of UO2 in 
the dry ore, g/cm"' 


7C/O2,i(ro) 


1.128 


1.680 


2.38 
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Present-day density of the dry 
rock with Pb, g/cm^ 




1.71 


1.61 


1.44 


6 


Density of H2O at T=300K 
(with water of crystallization) 


7h"2 


0.355 
0.455 


0.355 
0.455 


0.355 
0.455 
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Initial density of gj^U, g/cm^ 


78,«(0) 


1.305 


1.952 


2.758 
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Initial density of 92^U, g/cm" 


75,40) 


4.216-10"^ 


6.28-10"^ 


8.91-10~^ 


9 


Initial density of U, g/cm"' 


JuAO) 


1.347 


2.015 


2.847 


10 


Initial density of UO2, g/cm'' 


7(7O2,i(0) 


1.528 


2.286 


3.230 


11 


Initial density of ore, g/cm"* 


7.(0) 


3.51 


4.09 


4.78 


12 


Initial density of ore without 
UO2, g/cm" 


7o,«(0) 


1.98 


1.80 


1.55 


13 


Initial fraction of U in the 
dry ore, % 


YuAO) 


38.4 


49.3 


59.6 
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TABLE III: Specific weight ^k,i{0) and nuclear density Nk,i(0) in the compositions of three variants of the reactor core |44j 
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1.3683 
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1.3809 
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iH 


1 nn7Q 


3.97-10"^ 


2.373-10"^ 


3.97-10"^ 


2.373-10"^ 


3.97-10"^ 
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i4Si 


28.086 


0.398 


8.534- 10~^ 


0.362 


7.762-10"^ 


0.312 


6.690-10"^ 


7 




26.982 


0.270 


6.026-10"-' 


0.245 


5.468-10"^ 


0.211 


4.709-10"'' 
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12Mg 


24.305 


0.125 


3.097-10"^ 


0.113 


2.800-10"^ 


0.097 


2.403-10"^ 


9 


26Fe 


55.847 


0.286 


3.084-10"^ 


0.260 


2.804-10"^ 


0.224 


2.415-10"^ 


10 


19 K 


39.098 


0.0302 


4.652-10"* 


0.0275 


4.235-10"" 


0.0237 


3.653-10"" 






3.509 




4.081 




4.780 





of the cycle of the Oklo reactor. From eq.lO and eq.(|7I) we get for the ratio ju,i{Q)/lu,i{To) 

luAO)/luATo) = (1 - C5)exp(ro/T8) + C5exp(ro/T5) = 1.355, 



(8) 



independent of 7f7^i(To). The initial concentration Nk.i{0) of nuclei, which is needed for the calculations, was calculated 
from the formula 



NkAO)^lkAO)NA/Ak, 



(9) 



where Na = 6.022 • 10^^ mol"^ is the Avogadro number and Ak is the atomic weight. For ^gf U and ^gfU we have used 
the data of Tabled for the atoms of the rock we used the percentages by weight from Tabled The oxygen content of 
water was added to the oxygen of the core. The composition of the fresh core RZ2 that we used in calculations with 
different initial content of uranium is shown in Table ITTll Although the accuracy of the densities of some elements in 
this table is only a few percent, the values of Nk,i{0) are given with four decimal places for reproducibility of results. 

It follows from Table Cil that the enrichment of isotope ^sfu ((5(0) = 75(0)/7c/(0)) was (5(0) = 3.1% 1.8 billion 
years ago. Uranium of such enrichment is used in Russian VVER reactors of atomic power stations. Since the ratio 
of nuclei U/H is about equal and the sizes of both reactors are comparable, one can immediately and without any 
calculation say that a chain reaction was possible in Oklo jT?il |. 



Calculation of the fresh core 



1. Semianalytical calculation of core RZ2 



Consider first the bare reactor without reflector. Since the reactor is large compared with the neutron migration 
length M, one can apply the single-group diffusion theory [s^^^l- For a stationary neutron flux $(r) the following 
equation holds: 



1 

M2 



The solution that satisfies this equation with boundary conditions ea. pUf) is 

'2,405r' 



^(f) = $ocos( —x]Jo 



where Jo{Br) is the zeroth Bessel function, and the effective multiplication factor is 



1 + A'PB^ 



— -D// + i>/j , -Off 



TT 

H 



Br 



2.405 
R 



(10) 



(11) 



(12) 
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TABLE IV: Two-group parameters of the fresh cores of a cyhndrical bare reactor with different content of uranium i. Thickness 
of the core H = Im; average temperature in the core T = 300K; density of water 7H2O ~ 0.355 g/cm"^ 



Variant of the core 
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2 


3 


Relative initial density of U 














in 


the ore Yu,i{0), % 


38.4 


49.42 


59.6 


Radius of the active core, R, m 


9.9 


8.1 


6.8 


Computer code 


MCNP4C 


MCU REA 


MCNP4C 


MCU REA 


MCNP4C 


MCU REA 


1 




1.0965(1) 


1.0971(1) 


1.1238(2) 


1.1271(1) 


1.1247(2) 


1.1306(1) 


2 




1.1501(2) 


11/1 r\r\ /<~»\ 

1.1499(2) 


1.1750(2) 


1 1 /I \ 

1.1771(1) 


1 1 ^<~i 1 / 1 \ 

1.1721(1) 


1 1 '—7f'r\ / 1 \ 

1.1769(1) 
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Leakage of fast (/If), 
















(E > 0.625 eV) 


4.02T0"^ 


3.96-10"^ 


3.84T0~^ 


3.76T0~^ 


3.63(1)-10~^ 


3.55T0"^ 
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Leakage of thermal (Cth), 
















(E < 0.625 eV) 


0.434T0"^ 


0.43540"^ 


0.317T0"^ 


0.316T0"^ 


0.230(1)-10"^ 


0.230-10"^ 
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Square of diffusion length 
















= Ct,J[{l - Cth)Blo] cm" 


4.4 


4.4 


3.2 


3.2 


2.3 


2.3 


6 


Age r = - log(l - Cf)/Boo, cm^ 


41.3 


40.7 


39.3 


38.5 


37.0 


36.1 


7 


Total migration area 
















M'^ = K^T + L^, cm^ 


51.9 


50.2 


49.4 


48.5 


45.8 


44.8 


8 


B'io = tt''/H'' + 2.4048Vi^^ cm'^ 
TY^/H^ = 0.98696 • 10~^ cm"^ 


0.9929240"^ 


0.99584T0~^ 


0.99954-10"^ 
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K'^'f'f by Eq.Gl 


1.0937(2) 


1.0943(2) 


1.1199(2) 


1.1228(1) 


1.1208(2) 


1.1265(1) 
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Account of Eq. 11211 Apap — Peff, % 


-0.23(2) 


-0.23(2) 


-0.31(2) 


-0.34(2) 


-0.31(3) 


-0.32(1) 



In Table ITVI we show the two constants, i^oo.i and Mf, calculated with codes MCNP4C and MCU-REA for three 
different cores ^^lES- These constants are needed to calculate K^s,i by formula (|12|l . The values of -ftToo.i calculated 
for one and the same composition differ by a few tenth of a per cent; the values of M"^ — K^ iTi + differ by a 

few per cent. In row 9 of Table HVl we show the values of K^^^ calculated with the approximate formula H12|l . They 
are smaller than the direct calculations using Monte Carlo code (row 1). The difference in reactivity amounts to 
Api = —(0.2 — 0.3)%. The diffusion length in the fresh core is L = 1.6 — 2.1 cm, and the total migration length is 
Af = 6 — 7 cm. These lengths get less with increasing uranium concentration. 
The mean neutron flux, averaged over the reactor, is 

(Ji(2.405) — 0.51905). From formula H13|l we get the following formula of the volume nonuniformity coefficient Ky, 
independent of R and H: 

Kv = ^o/^^-- ^ = 3.638. (14) 

^ 2 2Ji(2.405) ^ ' 

This formula is useful to check the accuracy of calculation of the spatial distribution <I>(r). The absolute value of the 
mean neutron flux for Pp = 2.5 • 10"^ MW is equal to 

$ = (^-^Pp = 1.88- 10^5 n/s- (p, (15) 

where (p is the neutron flux per cm^ and one fast fission neutron which is calculated with the Monte Carlo code; 
Vf /Ef = 7.5 • 10^^ n/MW-s is the number of fast neutrons per second and a power of 1 MW {Ef is the fission energy; 
Vf \s the number of fast neutrons per fission). For thermal neutrons formula 115() holds with tp^^. The mean thermal 
neutron flux with energies En < 0.625 eV is very small in the case of Yu,2iO) — 49.4% it is $th = 0.63 • 10^ n/cm^-s. 
The thermal flux in the center of the core is <I'"jj'*'^(0) = 2.00 • 10^ n/cm^s. The total mean flux, integrated over all 
energies, is equal to $ = 3.9 • 10* n/cm^s. These results were found using code MCU-REA. The results of calculations 
using other Monte Carlo codes are similar (see Table lv|l . The low neutron flux determines the speciflcs of the function 
of the reactor. 

As a reflector one can assume the same core but without uranium. The analytical calculations for the reactor with 
reflector are more cumbersome. Therefore we used numerical methods for these calculations. The results are shown 
in Table IVll Both Monte Carlo programs give values of the reactivity reserve for variant 1 of the core which coincide 
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TABLE V; Calculations of the total flux ^tot averaged over the bare reactor, and of the average thermal flux 
^th {En < 0.625 eV) for the fresh core using three different codes. ^02(0) — ■i9A%U in the dry ore; i^^^o ~ 0.355; 
Pc = 0.1 MPa; r = 300K; Pp = 25 kW 





MCU-REA 


MCB 


MCNP4C 


K.ff 


1.1271(1) 


1.1237(1) 


1.1243(2) 


$iot, -10** n/cm"s 


3.9 


4.4 


4.4 


[En < 0.625 eV), -lO** n/cm^s 


0.63 


0.64 


0.63 



within the statistical accuracy. The difference of reactivity is 0.26% for core RZ2 of the bare reactor and 0.46% for 
core 3. For the reactor with reflector the difference is smaller: 0.20% and 0.39%, respectively. Since the migration 
length is small, a reflector of thickness A = 0.5 m is practically infinite: the results for a reflector of thickness 
A = 0.5 m coincide with those for A = 1 m. Compared with the bare reactor, the reflector makes a contribution of 
6p{A) = 0.8 0.9%. This contribution drops with increasing uranium content in the core. The cold reactor with a 
fresh core is strongly overcritical, since temperature and void effects have not yet been taken into account, also the 
initial strong absorbers which afterwards burn up rapidly. In Table IVlTl we show the neutron capture in the infinite 
fresh core per fast fission neutron. Capture by ^gfU amounts to 55.7% and by ^gfU to 33.8%. These are followed by 
hydrogen (3.9%), iron (3.8%), silicon (0.8%) etc. Code MCNP4C gives simflar values. 



TABLE VI: Calculation of K^ff and pi for cylindrical core of thickness Ho = Im and radius R for three different initial contents 
of uranium in the dry ore. Density of water 7H2O = 0.355g/cm''; T — 300K |43 



Variant of core 


1 


2 


3 


1 


Relative initial density of U in 
core YuAO), % 


38.4 


49.4 


59.6 


Speciflc weight of the dry ore 
7,(0), g/cm^ 


3.51 


4.08 


4.78 


Radius of active core R, m 


9.85 


8.07 


6.78 


Thickness of reflector A, m 


1 0.5 1 1 


1 0.5 1 1 


1 0.5 1 1 


2 


MCNP4C: 

Keff 

Reactivity pi, % 
<Sp,(A) =p,(A)-pi(0) 


1.0965(1) 
8.80(1) 



1.1077(1) 
9.72(1) 
0.92(1) 


1.1077(1) 
9.72(1) 
0.92(1) 


I. 1238(2) 

II. 02(2) 



I. 1357(2) 

II. 95(2) 
0.93(3) 


I. 1352(1) 

II. 91(1) 
0.89(2) 


I. 1247(2) 

II. 09(2) 



I. 1353(2) 

II. 92(2) 
0.83(3) 


I. 1357(1) 

II. 95(1) 
0.86(1) 


3 


MCU-REA: 

Keff 

Reactivity pi, % 
<Sp,(A) ^p^{A)-p,{0) 


1.0971(1) 
8.85(1) 



1.1080(1) 
9.75(1) 
0.90(1) 


1.1080(1) 
9.75(1) 
0.90(1) 


I. 1271(1) 

II. 28(1) 



1.1383(1) 
12.15(1) 
0.87(1) 


1.1381(1) 
12.13(1) 
0.85(1) 


I. 1306(1) 

II. 55(1) 



1.1407(1) 
12.33(1) 
0.78(1) 


1.1405(1) 
12.32(1) 
0.77(1) 



TABLE VII: Absorption < Ea4>V >fc and absorption with fission < S/^l^ >fc for infinite medium, normalized to one capture, 
for core variant Yu,i{Q) = 49A%U; T = 300K; % 





MCNP4C 


MCU REA 


< Ea<E>y > 


< E/$V > 


< Ea<E>V > 


< E/$\/ > 


1 




55.6 


45.9 


55.7 


46.0 


2 




34.0 


2.1 


33.8 


2.1 


3 





0.7 




0.7 




4 


H 


4.0 




3.9 




5 


Si 


0.8 




0.8 




6 


Al 


0.7 




0.7 




7 


Mg 


0.1 




0.1 




8 


Fe 


3.7 




3.8 




9 


K 


0.5 




0.5 




10 




100.0 


48.0 


100.0 


48.1 
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2. Power effect 



Reactor Oklo is controlled by the core temperature Tc ^46] . During heating the water was driven out of the core 
until the multiplication factor was equal to one. At first the large overcriticality was compensated by the power ejfect, 
which is the sum of the temperature and void effects. In Table IVlTll and Fig[21we show the dependence of the water 
density on the temperature for several pressures. At a pressure of 100 MPa in the Oklo reactor and Tc = TOOK, the 
density of water is 65% of its value for Tq — 300K and normal pressure. In this case the difference between crystalline 
and free water disappears apparently. The power effects is shown in FiglHl [13. Near Tc = 700K all the KcSi 
become equal to one. Therefore we can assume Tc — TOOK as the most likely temperature of the fresh active core 
(we neglect the small difference between the temperatures of fuel and water). For the variant of the composition of 
core RZ2 , the power effect is App — —11.6%. The void effect at 700 K accounts for 73% of this value, and the 
temperature effect for 27%. These results were obtained with code MCNP4C. Code MCU-REA gives similar values. 
In Table |IV| we show the numerical values of Kcs{Tc), Koo{Tc) and AP{Tc), calculated with code MCU-REA for 
the bare reactor with core composition i — 2. With increased temperature Koo{Tc) drops and and the leakage 
increase (Fig^J. 



TABLE VIII: Temperature dependence of the water density at a pressure of Pc ~ 100 MPa in the core 



T{K) 


IH2O, g/cm^ 


1^17^20, g/cm'^ 


i^2'yH20, g/cm^ 


300 


1.037 


0.368 


0.472 


350 


1.014 


0.360 


0.461 


390 


0.994 


0.353 


0.452 


400 


0.982 


0.348 


0.447 


500 


0.900 


0.319 


0.409 


536 


0.864 


0.307 


0.393 


540 


0.859 


0.306 


0.391 


600 


0.792 


0.281 


0.360 


700 


0.651 


0.231 


0.296 


800 


0.482 


0.171 


0.219 


900 


0.343 


0.122 


0.156 



) = 0,355; uj2 = 0,455. 




FIG. 2: Dependence of the water density 7H2O on the temperature T for different pressures P. 



To determine the spread of results depending on the uncertainty of the initial composition of the core the calculations 
were carried out over a wide range of the content of uranium {Yu,i{0) = 39.4 — 59.6% by weight) and of water 
{uj^u Q — 0.355 — 0.455) in the ore. For reliability the calculations for the bare reactor were carried out with two 




FIG. 3: Power effect 0, Dependence of the reactivity p (in %) of the fresh core of the Oklo reactor on the temperature 
Tc at a pressure of Pc = 100 MPa for three different initial compositions of the active core and different proportions of water: 
1 - 49.4 voL% U in ore, uj°„^o = 0.455; 2 - 49.4 vol.% U in ore, J\i^o = 0.405; 3 - 38.4 voL% U in ore, w^^o = 0.355. The 
calculations were done using code MCU-REA. 



codes: MCU-REA and MCNP4C. The results were additionally controlled by the single-group formula H12|) with the 
parameters shown in Fig^l The calculations of the Tc dependence of i^Toff for variants 2 and 3 are similar (Fig|2I). 
For variant 1 with lower uranium content (F[/i(0) = 38.4% by weight) the curve K^'g^c) lies visibly lower. For a 
water content of ^^20 ~ 0.455 the curve for variant 2 lies significantly higher. As a result the core temperature at 
which the reactor became critical was Tc — 725K with a spread of ±55K. Taking account of the fuel burn-up and of 
slogging of the reactor a loss of reactivity takes place. This leads to a drop of Tc- The calculations of burn- up are 
continued at present. 



,20 
18 
16 
14 
12 
10 
,08 
1,06 



00^ C ' 



_ 9 




• 

Okl 


- 49.4% U 




3 




a 


— 

6>^^^ — □ 


— K 

- 


^^^^ 








..■ g 
















B 








1 










^[ 

■j-^^a^''^ 






!^^^^^^ 




a 


3 — ' 






9 


3 



75 



70 



65 



60 



55 



50 



45 



300 



400 



500 



600 



700 



FIG. 4: Tc dependence of K^{Tc) and M'^{Tc) for a core containing 49.4 volume % of uranium in the present-day dry ore 
and a water content of uj%^o = 0.355 at Tc = 300K and Pc = 0.1 MPa. Calculations using code MCU-REA. ,44. .45.1 
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The reactor could have worked also in a pulsating mode: when the temperature exceeded 710 K, then the unbound 
water was boiled away and the reactor stopped on account of the void effect. Then the water returned and it started 
to work again 0, 0, ^ . However a detailed analysis of the pulsating mode of operation of the reactor is outside 
the scope of the present paper. 

3. Computational problems in calculations of large reactors 

When making Monte Carlo calculations of the neutron flux in large reactors, one encounters certain difficulties. In 
such reactors many generations are produced before a neutron that was created in the center of the reactor reaches 
its boundary. This time depends on the relation between the migration length and the size of the reactor pIs] . At 
Tc — 300K these values are for core i = 2 equal to M = 7 cm and R = 8.1 m. 230 generations are needed before 
a centrally produced neutron reaches the boundary, detects the boundary condition and returns to the center. In 
order to reproduce the spatial distribution of the neutron flux with sufScient accuracy one must calculate tens of 
such journeys. Analyzing the solution of the time dependent diffusion equation one finds that over 6000 cycles are 
needed to get the fundamental harmonic with an accuracy of a few per cent. Experience with such calculations shows 
that one needs (5 — 10) • 10'^ histories per cycle in order to keep an acceptable statistical accuracy. Thus we needed 
(4 — 6) • 10^ neutron trajectories for our calculations. For several hundred calculations we have explored an order of 
10^*^ trajectories taking up several months of continuous work of a modern PC cluster. In spite of such a large volume 
of calculations we could not find the volume nonuniformity coefficient Ky of the neutron fiux with good accuracy. To 
do these calculations we had to divide the core into tens of volume elements which led to a reduction of the statistical 
accuracy in each of them. As a result the value of Kv in formula (|14|l was reproducible with an accuracy not better 
than 10%. This is obviously insufficient to carry out the calculation of the burn- up that depends on the magnitude 
of the absolute ffux in different parts of the core. One must admit that the Monte Carlo method is not suitable for 
the calculation of large reactors and one must resort to different approaches. 

The reactor neutron spectrum below 0.625 eV is needed in order to average the cross sections of strong absorbers 
(e.g. ^Q2^m). The spectrum for three compositions of the fresh core without reflector, calculated with code MCNP4C 
for Tc — BOOK, is shown in Fig[Sl One can see small peaks which correspond to excitations of rotational and 
vibrational levels of II2O. For comparison we also show in Fig. [Slthc Maxwell neutron spectrum that was used by all 
previous authors to average the ^|fSm cross section 7, 8, 9, 10, 11, 12, 13]. The spectra are significantly different. 
The Maxwell spectrum has a much higher peak but is exponentially small above 0.3 eV where the reactor spectrum is 
a Fermi distribution. In our calculations we have used the Nelkin model of water which automatically takes account 
of the chemical bond of hydrogen nuclei. Calculations at other values of Tc yield similar results. 

III. VARIATION OF FUNDAMENTAL CONSTANTS 

The Oklo reactor is an instrument that is sensitive to the neutron cross sections in the distant past. Comparing 
them with current values one can estimate how constant they, and hence also the fundamental constants, are in time 

mum 

A. Early approaches 

In 1935 E. Miln posed the question: how do we know that the fundamental constants are actually constant in 
time '5^. He thought that the answer could be found only by experiment. A little later D. Dirac proposed that 
originally all constants were of one order of magnitude but that the gravitational constant dropped at a rate of 
G/G ~ —tg^ during the lifetime to of the univer se Isil Is^ . In 1967 G. Gamov suggested that, on the contrary, the 
electromagnetic constant is increasing: a/a ^ to |53l |. Both hypothesis were wrong since they contradicted geological 
and paleobotanical data from the early history of the Earth. Without entering into a detailed discussion of these and 
many other later publications on this subject one must admit that there is a problem of the experimental limit on 
the rate of change of the fundamental constants (see the early review by F. Dyson [s^l. 

The authors of Refs. 0, H, 0, noticed that the sensitivity to variations of the nuclear potential increases by 
several orders of magnitude if one considers neutron capture. Owing to the sharp resonances of the absorption cross 
section the nucleus is a finely tuned neutron receiver. A resonance shifts on the energy scale with changes of the nuclear 
potential similarly as the frequency of an ordinary radio receiver shifts when the parameters of the resonance circuit 
are changed (Fig|Hl) 55]. Qualitatively one can understand the absence of a significant shift of the near-threshold 
resonances on the grounds that all strong absorbers are highly burnt up in the Oklo reactor and weak absorbers are 
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700 



E^, meV 



FIG. 5: Neutron spectrum in the bare fresh core RZ2 at different initial uranium concentrations (water content u!%^o ~ 0.355). 
Calculations using code MCNP4C |i3,|i3. 

burnt up weakly (Fig[7|) 13, Is^. Holes in the distributions are seen for strong absorbers: ^Ij^m, ^gsEu, ^gl^d, ^gl^*^- 
The depth of burn-up, calculated using the present absorption values, are in satisfactory agreement with experiment, 
particularly if one remembers that the neutron spectrum over which one must average the cross section is not known 
very well. Thus in the 1.8 billion years since the work of the Oklo reactor, the resonances (or, in other words, the 
levels of the compound nuclei) have shifted by less than AEr ~ r7/2 {T^ = 0.1 eV). Therefore the average rate of 
the shift did not exceed 3 • 10~^^ eV/year. This value is at least three orders of magnitude less than the experimental 
limit on the rate of change of the transition energy in the decay of ^^"^Re |53 |. 



At present there are no theoretical calculations giving a reliable connection between the positions of all resonances 
with parameters of the nuclear potential. But already the preliminary qualitative estimates allow one to reduce the 

limits on the rates of change of the coupling constants of the strong and electromagnetic interactions a/ a and Sa/a. 
We confirm the absence of a power or logarithmic dependence on the lifetime of the universe. It is desirable to have a 
more detailed calculation of the influence of variations of the fundamental constants on the parameters of the neutron 
resonances. 



B. Basic formulae 

1. Averaging the Breit-Wigner formula 

When a slow neutron is captured by a nucleus of isotope ^gfSm, then a nuclear reaction takes place with formation 
of an excited intermediate compound nucleus and subsequent emission of m 7-quanta: 

n + i^^Sm ^ ^I^Sm* ^ ^I^Sm + m • 7 . (16) 

Near a strong S resonance one can neglect the effect of the other resonances and describe the cross section with the 
Breit-Wigner formula 
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B„=6-8 MeV 



/////////// 



n+ Z ^ Z 




FIG. 6: Strong absorber as a sensitive detector of a variation of Er |53 |. Shown on the left is the energy level density of the 
compound nucleus n Z ^^'^^ Z* . Shown on the right are resonances in the cross section of the reaction . The capture 
cross section behaves like ~ (F^/Er)^, where Er is the distance from the resonance and Vj is its width. Neutron capture is 
strongly affected by a shift of AEr- 



10" 



10' 



10' 



10-' 



10^- 



loS 



10- 



— I — 
144 



Eu 



"Sm 



o 

X 



X - Theory 
o - Experiment 



"Gd 



— I — 
156 



148 



152 



160 



FIG. 7: Comparison of the calculated (crosses) and measured( circles) concentrations of fission fragments Na relative to the 
'^^^Nd content for one of the samples of the Oklo reactor [loll5q|. 

where go = (2J + 1)/{2S + 1)(2/ + 1) is the statistical factor, S* = 1/2 is the electron spin, / is the nuclear spin and 
J is the spin of the compound nucleus. The full width is Ttot — r„(£') + F^, where Tn{E) and T-y are the neutron 
and 7 width, respectively. The neutron width is given by 57] 



V E() 

The parameters of the lowest resonances of a number of absorbers is given in Table ^ 



(18) 



In formula H17|) the neutron energy is given in the cm. frame: Ec — \mn \ Vl — Vk p. It depends on the velocities 
of the nucleus Vk and the neutron Vl in the lab frame and on the reduced mass m„. The reaction rate NkUj^kiEc) ■ Vc 
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TABLE IX: Lowest resonance parameters of strong absorbers \5T\ 





62 sm 




84' Gd 


Resonance energy (E,.), mcV 


97.3 


26.8 


31.4 


Neutron width at 1 eV (r"),incV 


1.71(3) 


0.635(10) 


2.66(5) 


Gamma width (F^), meV 


60.5(6) 


108(1) 


106(1) 


g - factor 


9/16 


5/8 


5/8 


Doppler width (Ad), meV at T = 700 K 


12.5 


6.4 


6.9 



with cross section H17(l and for an absorber of nuclear density Nk must be averaged over the nuclear spectra fk{Ek) 
and the neutron spectrum n{EL) (all spectra are normalized to one). The inverse nuclear burn-up time in an arbitrary 
point of the core is given by 



A^,fc(r) = A^fc dpkdpLfk{Ek)n{EL)a^^k{Ec)Vc ■ 



(19) 



At high temperatures the gas approximation is valid for heavy nuclei of the absorber. Changing to integration over 
the cm. energy Ec and the neutron energy El and assuming a Maxwell nuclear spectrum, we get 



A^,fc(r) = Nk J n{EL)a^^k{Ec)VcF{Ec ^ EL)dELdEc , 



(20) 



F{Ec — ^ El) is the transformation function from the cm. to the lab system (for details see Ref. jSq): 



F{Ec ^ El) 



{A + l) 



2,/ttATEl 



exp 



A 
T 



1 



il + j)-Ec 



exp 



A 
T 



1 



'Er 



(21) 

and A = MA/TTin is the mass of nucleus A in units of the neutron mass. 

Close to a resonance we can neglect the second term in Ea. (|21|l and evaluate the first term in integral H2U|) by the 
saddle-point method. As a result the integral (|20|l takes on the following form in the vicinity of a resonance: 



A^,fc(T) = iVfe| (l + ^ jdEca^^k{Ec)Vc J dELn{EL)V 



El 



where the Gaussian 



r 



El-[1 + \]Ec 



is normalized to one and the Doppler width is equal to 

An = 



exp 



[EL-{l + ^)Ecy 



'^ElT 


1/2 


-AEcT 


1/2 


A 




U+iJ 





(22) 



(23) 



(24) 



The values of the Doppler widths for T = 700K are shown in Table Hxl Since all A^j ^ F^, function H23|l can be 
replaced hy 5 {\ El ~ AEcl{A + 1) |) and integral becomes 



A 



o--t.kiEc)Vcn 



A 



Er 



dE, 



c ■ 



(25) 



The correction 2/ A is of magnitude 1%. If the neutron spectrum is Maxwellian in the cm. frame, then it is also 
Maxwellian (with reduced neutron mass) when the nuclear motion is taken into account. It can be shown that Ea. H25() 
is valid at all energies if the distribution of nuclei and neutrons is Maxwellian js^ . 

Therefore it is not surprising that the authors of Ref. did not notice any deviations from formula (|25l) in their 
numerical calculation that took account of the thermal motion of the target nuclei (Doppler effect). However, the 
situation is different if the neutron spectrum is not Maxwellian. In this case one must use formula H22|l instead of the 
simple formula H25|l . 
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To average the capture cross section of samarium one normahzes the cross section, integrated over the neutron flux 
spectrum {n{E,T)v), traditionally not by the integrated flux but by the product of the velocity — 2200 m/s and 
the integrated neutron density {n{E)) [Ta lT^l5fll |: 



o-7,fc(r) 



/ (y-,,k[EL)n{EL) ■ VbdEL 
vOJn{EL)dEL 



(26) 



If the cross section a^^kiEL) has a. 1/vl behaviour, then the integral of a^^kiT) is constant. From formula one 
has 



0'7,fc(r) 



4T 
ttTo 



(27) 



where Tn = 300 K = 25.9 meV. Useful is also the relation lie 



(7$ — (7$, where 



2 T 



(28) 



We have evaluated the cross section o'^^SmiT) of ^|2Sm without recourse to any approximations. For 



JdEkdELfs,n{Ek) (Ec) ■ Vcn{EL) 
v°JdELn{EL) 



(29) 



For the calculations we used the computer package MATHEMATICA In FigEI 

we show the values of 

o'-y,Sm{T, AEr) at six temperatures T = 300 — lOOOK for a shift of the resonance position AEr = ±0.2 eV. The 
curves have a maximum at negative shifts of the resonance; the maximum of the curves is higher at lower temper- 
ature T. At the point AE^ ~ and at T = 293K the cross section calculated as the contribution of the closest 
resonance is equal to cr^.gm (293 K) = 39.2 kb. The contribution of higher positive resonances is (T^g„j(293K) = 0.6 
kb and negative ones is cr^^^ (293 K) = 0.3 kb ^7}. The total cross section (as measured on a neutron beam) is 

tT*°^^ (293 K) = 40.1 kb. At small energy shifts AEr (^^,sm ^nd cr:^^g„^ practically do not change. In Refs.p H ll H l^ 
the total cross section 40.1 kb has been used instead of the single resonance one 39.2 kb. Therefore the curves in Refs. 
in El m are higher by 40.1 kb/39.2 kb, i.e. by 2.5%. 
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FIG. 8: Dependence of the cross section a^sm, averaged over a Maxwell neutron spectrum, on the resonance shift AEr and 
on the temperature: T — (300 — 800)K ea. ll29ll . The curves are for fixed temperatures with intervals of 100 K. The upper left 
curve is for T = 300K. 
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2. Taking account of the reactor spectrum 

In Fig.?? we show the resuhs of calculating a^,sm{Tc t ^Er) with the Maxwell spectrum replaced by the reactor 
spectrum nfi{E]^^Tc)- The central curve 2 is the result of the calculation using code MCNP4C for the fresh core with 
^(72(0) = 49.4%U in the dry ore and with w^^o ~ 0.405 at T = 725K. For comparison we also show the cross section 
averaged over the Maxwell spectrum at T = 725K for the same composition of the core (curve 4)- Curves 4 and 2 are 
significantly different, especially at negative AiJ^: curve 2 lies distinctly lower. The maximum of curve 2 is 1.5 times 
lower than the maximum of curve 4- At lower temperatures this difference is even greater. Thus we conclude that we 
have proved a significant effect of the reactor spectrum on the cross section of ^||Sm. 
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FIG. 9: Dependence of the thermal neutron capture cross section of isotope 62 Sm on the core temperature Tc and on the 
resonance shift AEr'. (T^,sm{T, AEr) |4^ |45|. The curves are for cross sections, averaged over the reactor spectrum of the 
fresh core at three different initial states: 1 - Yu2{Q) = 38.4 vol.%U, wh^o = 0.355, Tc = 670K; 2 - Yu2{0) = 49.4 vol.%U, 
UJH20 = 0.405, Tc = 725K; 3 - Yu2{Q) = 49.4 vol.%U, ujh^o = 0.455, Tc = 780K, Pc = 100 MPa. For comparison we show 
the cross section averaged over the Maxwell spectrum (curve 4) for initial composition Yu2{0) = 38.4 vol.%U, UJH2O = 0.405, 
Tc = 725K. Shown is the error corridor of measured values d'^'^g^ fTlLll2|. 



In order to determine the dependence of a^^sm{Tc, A.Er) on the uncertainty in the initial active core composition, 
we have calculated the values for the two outermost curves of Fig|51 Curve 3 of this figure corresponds to an initial 
content of Yiji{0) = 38.4%U in the ore, w^^q = 0.355 and Tc — 670K; curve 1 corresponds to an initial content of 
49.4%, UJH2O ~ 0-455 at T = 780K. Since the numerical constants are known only for values of Tc which are multiples 
of 100, we have done the calculations for Tc = (600, 700, 800)K and interpolated to intermediate temperatures. The 
broadening of curve 2 on account of the scatter of temperatures is small. Experimental data of d'^^g'^{T) for core 2 
are presented in Ref. (Table Ixjl (see also Refs. [Q). The labeling of sample SC36-1418 indicates that the sample 
was taken from bore-hole SC36 at a depth of 14 m 18 cm. The mean value a^ g,^ — (73.2 ± 9.4) kb is shown in FigEl 
Curve 1 (Tc = 670K) crosses the lower limit of a^^^sm — 64 kb to the left of point AiJ^^' = —73 meV, and curve 3 

(2) 

{Tc = 780K) to the right at AE}- ' = +62 meV. The possible shift of the resonance is therefore given by these limits: 



-73meV < AEr < 62meV. 



(30) 



3. Connection between AEr and Sa/a 



The shift AEr must be related to a variation of the fundamental constants, for instance to a shift of the electro- 
magnetic constant a = 1/137.036. This has been done by Damour and Dyson The change of the Coulomb 
energy contribution AHc to the energy of the level in the nuclear potential AEc, that results from a change of a, is 
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TABLE X; Experimental values of asm for 15 samples from Okie reactor core RZ2 







(5-Sm, kb 




1 


KN50-3548 


93 


I56| 


2 


SC36-1440 


73 


I61| 


3 


SC36-1410/3 


73 


I61| 


4 


SC36-1413/3 


83 


[61] 


5 


SC36-1418 


64 


I61| 


6 


SC39-1383 


66 


[27, 61] 


7 


SC39-1385 


69 


[27, 6ij 


8 


SC39-1389 


64 


[27, 61j 


9 


SC39-1390 


82 


[27, 61j 


iU 




QO 
OZ 


127, (Jil 


11 


SC39-1393 


68 


[27, ^ 


12 


SC35bis-2126 


57 


[27, 6JJ 


13 


SC35bis-2130 


81 


[27, 61J 


14 


SC35bis-2134 


71 


[27, 61] 


15 


SC52-1472 


72 


[32] 




(TSm ± AlTSni, kb 


(73.2±9.4) 





given by 



^AEc^liAHc) , (31) 



In first perturbativc approximation the dominant contribution to AEc is given by the isotopic effect {^^, p. 568) 

AEc = En^Y\ ^<^(0) \^ Z-e^<R^> . (32) 

Here "06(0) is the wave function of the s wave electrons in the nucleons and (i?^) = {Z ■ e)^^ J pR^dV, where p{R) 
describes the proton charge distribution in the nucleus. Damour and Dyson have estimated the value of (i?^) for the 
excited nucleus ^Ij^m* from the neighbouring isotopes. They found 

7W = ASfl = -(l.l±0.1)MeV. (33) 

Combining this value with the shift of the resonance AEr in formula H3U|) . we get for (3AEr/Ai 

-5.6 • 10"^ < ^/a < 6.6 • 10"^. (34) 

Because of the negative value of A4 the limits on 6a/ a change their places. For the past time (—To) the product 

(— TqA^) is positive and hence the limits on Sa/a = A£',.(— TqAI) are restored to their previous places. Note that 
traditionally Sa/a is defined by Sa = {aokio — ctnow) / ct- This shift of a lies in a narrower range than in Ref. |llL Il2| . 
Assuming a linear change of the e.m. constant during the time To, we get the following limit on the relative rate of 
change: 

-3.7 • 10-" year-i < fc/a < 3.1 • 10"" year'^ . (35) 

Thus, within the limits given by Eq. H35|l . the e.m. constant changes for the fresh reactor Oklo core with zero speed, 
i.e. it remains constant. 



C. Review of previous work 

The work of Shlyakhter (1976, 1983) SHE/- The authors of Refs. 0, S were the first to point out the 

possibility of using the data of the natural nuclear reactor Oklo to find the most precise limits on the rate of change 
of the hmdamental constants. The most convenient data are those of the strong absorbers, e.g. of ^g2Sm. For this 
isotope Shlyakhter calculated at T = 300K the dependence of the change of the cross section on the resonance by an 
amount of A£:^: a-y^sm{Tc, AEr) (FigUHl) 
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FIG. 10: Change of the thermal (T — 25meV= 273K) capture cross section of ^62 Sm for a uniform shift of all resonances by 



He compared this curve with the experimental data available at the time: cr^^sm ~ (^^ ^ ^) ^ possible shift of 
the first resonance within two standard deviations (95% confidence level) was found to be 

6Ef^^ < 20 meY . (36) 

(In going from a.y^sm{T, AEr) to ay,sm{T, AEr), all values must be multiplied by 1.18 according to Eq. H27|l and 

^^^Sm = (65 ± 9.5) kb, but this does not affect the value of 6Ef2'^). To estimate A4, Shlyakhter used data on the 
compressibility of the nucleus; he found A4 — —2 MeV 0- With a linear dependence of the change of a with time 
for To = 2 • 10^ years the limit on the rate of change of a is 

Sa/a < 0.5 • 10"^^ year^^ . (37) 
Using the present, more accurate value of ~ — 1.1 MeV, we get 

^/a < 1 • lO^^Vear-^ (38) 

It should emphasized again that this limit was found for only one temperature: T = BOOK. 

The work of Petrov (1977) f73/ . In this paper the resonance shift 5Er was estimated from the shifts of the widths of 
a few strong absorbers. The resonances of strong absorbers lie close to a zero energy of the neutron, and the resonance 
energy is of the capture width: 5Er ~ ~ 0.1. The capture cross section of these absorbers for thermal neutrons 
changes sharply when the resonance is shifted by an amount of the order of The analysis of experimental data 

for Sm and ^IsEu, taking account of a threefold standard deviation and the uncertainly of the core temperature, 
shows that the shift 5Er of the resonance since the activity of the Oklo reactor does not exceed ±0.05 eV [l^. The 
results of the measurement of the concentration of rare earth elements with respect to "'^go-'^'i {the second branch of 
the mass distribution of the fission fragments) in one of the Oklo samples are shown in Fig[7| A more conservative 
estimate in Ref. lOj is | 5Er |< 50 meV, i.e. 2.5 times higher than Shlyakhter's estimate. Using the modern value 
M ~ -1.1 MeV, we find for fo/a < '5E'r/\ M \ 

fc/a =Ie;/| M |< 4.5- 10~^ (39) 
This is almost 5 times greater than Shlyakhter's optimistic estimate. For the rate of change Sa/a we get 

la/a < 2.5 • 10"^^ year"! . (40) 

This is less by a factor of 2 than the limit found 20 years later by Damour and Dyson |^ 0| . The reason of this 
discrepancy is the use of only one temperature {Tq = BOOK). Although the dependence of AEr on Tc was noted in 
Ref. no calculations of the effect of the temperature were carried out. 
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The work of Damour and Dyson (DD) (1996) 0, CI/- The dependence l^Er{Tc) was analysed 20 years later in 
the paper DD pHIT^ . They have repeated the analysis of Shlyakhter and came to the conclusion that it was correct. 
DD also updated Shlyakhter's data in three directions: 
(i) They employed a large amount of experimental data (see Table Ixj) . 

{ii) They have taken account of the great uncertainty of the reactor temperature, Tc = (450— 1000)°C (Fig. Illfl . As 
a result they made a conservative estimate of the mean shift of the resonance 



-120 meV < /^Er < 90meV. 

The range of the shift AEri — AEr2 = 210 meV is 1.5 times greater than the range of the shift in our paper. 
(in) They have calculated the value M = -(1.1 ± 0.1) MeV, but used M = -1 MeV. For Sa/a DD found 

-9.0 • 10"* < < 12 • 10"* . 

This leads to the following limits on the rate of change Sa/a: 

-6.7 • 10"" year"i < S^/a < 5.0 • 10"^^ year"^ . 



(41) 



(42) 



(43) 



Since ^If^m burns up 100 times faster than ^gfU, therefore the only ^If^m found in the stopped reactor is that which 
was produced immediately before the end of the cycle. As a consequence DD emphasized that one must know the 
detailed distribution of the nuclear reaction products of the end of the cycle to make a detailed analysis. 
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FIG. 11: Dependence of the thermal neutron capture cross section on the resonance shift AEr for the two temperatures of the 
Maxwell neutron spectrum {Damour & Dyson [lllll^'). The energies of the crossing of the lower limit of data are AE^[ — —120 
meV and AE^i = 90 meV. 



The work of Fujii et al. (2000, 2002) fJl. fg^/ . The experimental data on the measurement of a^^k of strong 
absorbers were presented in the papers of Fujii et al. Of five experimental points, four are from core RZIO (Table 
IX1() and one from another core, RZ13 , and therefore we have omitted it. Core RZIO lies at a depth of about 150 m 
from the surface of the quarry. As we do not have any detailed data on the size and composition of core RZIO , we 
shall assume them to be similar to those of core RZ2 . For ^IfSm the mean value of the four points of Table IXj is 

= (90.7 ± 8.2) kb. (44) 

This value is noticeably greater than for core RZ2 [(72.3±9.4) kb (see subsection llll B Q^ . The error bars of both values 
do not even touch and their difference remains significant. For increased reliability the number of measurements for 
core RZIO should be increased. Possibly this core has finished its cycle at a lower temperature. 

In FigElwe show the dependence oia-^^sm{Tc , AEr) on the shift of the resonance for a Maxwell distribution 0,1^. 
The authors have estimated (on the grounds of indirect considerations) the uncertainty in Tc — (180 — 400)°C = 
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TABLE XI: Experimental values of rj^^^ (kb) for strong absorbers, measured in the middle of the core RZIO , located at a 
depth of 150 m 









54' Gd 






afS, kb 


kb 


<7^57, kb 


1 


SF84-1469 


83.6 


30.9 


83.3 


2 


SF84-1480 


96.5 


16.8 


8.0 


3 


SF84-1485 


83.8 


17.8 


14.3 


4 


SF84-1492 


99.0 


36.7 


73.7 


5 


CTfc 


90.7 


25.6 


44.8 


6 




8.2 


9.8 


39.2 



(453 — 673)K. They also show the experimental data of formula 1)44(1 . The intersection of the limiting curves with the 
lower limit cr^^g^ = 82.5 kb yield the following possible shift of AiS^: 



-105 mcV < A£;,. < +20meV. 



(45) 



In Fig we also show for comparison the curves a^^sm{Tc: ^Er) for the reactor spectrum of the fresh core at 
Tc = 400°C, Yui{0) = 38.4%U in the ore, u^^^ = 0.355 and P = 100 MPa. The possible shift of AS^ for the 



experimental data of formula (|44|) lie in a narrower interval than in the paper DD 

-120 meV < l^Er < 20meV 
From Eq.(@ni and using M = —1.1 MeV we get 



Ltie pos 

Urn 



-1.8 ■ 10"* < Sa/a < 9.5 • 10" 



and 



-5.3 • 10"^^ year^i < Sa/a < 1.0- 10"" year~^ . 
Thus in this case too we do not find with certainty a nonzero deviation of the change of a. 
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(46) 



(47) 



(48) 



IsE , meV 

FIG. 12; Dependence <TT.,Sm(Tb, AS^) for the Maxwell distribution at Tc = (450, 570, 670)K {curves 1-3) (Fujii et al.[i3,E3) . 
The data correspond to Table El ct^^^ = (90.7 ± 8.2) kb, AE^i = 20 meV. For comparison we show curve 4 for the reactor 
spectrum at Tc = 670K; F!7i(0) = 38.4 vol.%U; ^loh^o = 0.355; Pc = 100 MPa. 

The work of Lamoreaux and Torgerson (LT) (2004) El/- These authors noted (two years after Ref. [l^ ) 
that the reactor spectrum contains in addition to the Maxwell tail also the spectrum of the moderated neutrons 
(Fermi spectrum). They averaged a~f^smiTc, AEr) at Tc = 600K over this spectrum and, comparing this curve 
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with the experimental data of Fujii et al. [b^, found a shift AE'r = (— 45^]^5) meV. Thus they found a shift of the 
cross section of ^gaSm. Let us consider the LT model in more detail in order to understand this result. The age 
of the Oklo reactor is Tq 2 • 10^ year. The relative content of ^sfu is ^5(0) = 3.7% (and not 3.1%). The ratio 
of H to U is fn = Nh/Njj = 3. The cross section of the burning up admixtures (e.g. lithium) per atom of U is 
Pu = ^i'^a/^u = 2 b. The ratio of the thermal neutron capture cross section to the slowing down cross section 
of epithermal neutrons is A = J^aiVr / {^s fiA) = 2. The temperature of the core is Tc = 600K= 327°C. 

The following comments are appropriate concerning these parameters. The age of the reactor is 1.8 • 10^ year and 
not 2 • 10^ year. The value of ^5(0) = 3.7% yields cross sections £5,0(0) and T,5j{0) which are too large by 1.2 times. 
The ratio Jh = Nh/Nu = 3 holds approximately for YusiO) = 59.6% and w^^^ = 0-355 (see Table Hflt. For the 
condition fn = 3 to hold exactly one must reduce ^^^20 ~ 0-355 to w^^o ~ 0.323. The concentration of uranium 
nuclei can be kept in the calculations at its former value: Njj = 0.7205 • 10~^ U/cm-b. The absorption cross section 
of 3Li nuclei per uranium nucleus, (3 = NLi(Ja.Li/Nu = 2 b at T = 300K and normal pressure results in the following 
concentration of lithium nuclei: N^i = 2.088 - lO""* Li/cm-b. 

The ratio of the capture cross section of thermal neutrons to the scattering cross section of epithermal neutrons is 
too large in the LT paper. At T = 300K, LT use for the capture cross section the value (Ja,kNk/Nu = 31.1 b/U. 
The scattering cross section of a free hydrogen nucleus is 20.5 b "sT] (for bound hydrogen it is greater). At = 3 
the value of A is A = 2 • 31.1/3 • 20.5 — 1.01 and not 2. Even though, in repeating the calculation of LT we have used 
the value A = 2. The curve CF^^sm{Tc, ^Er) is shown in Fig. US (curve 1). The value of (J^J„(r) = (90.7 ± 8.2) kb 
are taken from Table IXll Recall that Table ITTl contains only four experimental points instead of five. This results in 
a change of cr^^g^^ and its error as compared to LT . Curve 1 was obtained by interpolation of curves 2 for A = 1 

and for A = 2 to the value A = 2^300K/600K = ^/2. We have found a negative value for the energy Ai?r2 = — 24 
meV, closer to zero than the value A£'r2 = —38 meV in LT . From Fig^jone can clearly see a specific feature of the 
result of LT . It is enough to take A > 1.41 for the curve not to intersect the error corridor, and for A < 1.41 the 
shift Ai?r2 is strongly reduced. 
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FIG. 13: Dependence a-jSm{T, AEr) at Tc ~ 600K for the neutron flux of Lamoreaux & Torgergerson Ref. |64l| . calculated by 
us for A = \/2 ( curve 1). For comparison we show results for A = 1 {curve 2) and A = 2 (curve 3). Curves 4 and 5 are our 
calculations for the core composition of LT; curve 4 is without the power effect (PE); curve 5 with PE (computation with code 
MCU-REA 



At a stiffness parameter A ^ 1 the spectrum is distorted on account of large absorption by the strong absorbers at 
small energies. Under these conditions the spectrum cannot be considered to be Maxwellian. It can be found only by 
direct Monte Carlo calculation. We have carried out the calculation of the reactor spectrum and of the distribution 
for a composition of the core as described in items 3-5, using the code MCU REA. In the LT paper no absolute value 
of Nij{0) was given, and we choose Njj = 0.7205 • 10^^ U/cm-b. The calculations were done both without taking 
account of the power effect {PE) (curve 4) and with taking account of the PE (curve 5). Both calculations cross the 

error corridor at l^E^^ > {^E^^ = 4 meV). This means that in this case Sa/a = as well. 

In the LT model the neutron balance is maintained on account of a compensation of the fuel burn-up by the burn-up 
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of strong absorbers. This is far from reality. Since strong absorbers burn up faster than ^gf U, such a balance exists 
only at the beginning of the cycle. At the end of the cycle no strong absorber is left. The fast burning up strong 
absorber -'^^^Sm that is still present to this day was formed only at the end of the cycle when the LT model does not 
work any more. 



The results on a possible change of a based on the analysis of the cross section of 



62^ 



the Oklo reactor 



are summarized in Table IXIII For comparison we have included the cosmological results fFig ll4ll [b^ l68j | and the 

results of laboratory measurements |70| . A review of a possible change of the fundamental constants (experiment and 
theoretical interpretation) was recently published by Uzan 71, 72]. All results show that there are no grounds for an 
assertion that the e.m. constant has changed in the distant past. However there is a possibility that this conclusion 
will be revised when the fuel burn- up is taken into account. 



s 




FIG. 14: The data of Chand, Srianand et al. (filled circles) are plotted against the redshift f6^. Each point is the best 
fitted value obtained for individual systems using minimization. The open cirle is measurement from the Oklo reactor. 
The weighted average and la range measured by Murphy et al. [69l | are shown with horizontal dashed lines. Most of Chand 
et al. measurements are inconsistent with this range. The shaded region represents the weighted average and its 3a error: 
< 5a/oL >w= (-60 ± 60) ■ 10"**. Within 3a there is no variation of fundamental constants within the limits: 
-25 ■ 10"^^ year"^ < (Aa/aAt) < 12 • 10"^^ year"^ 



IV. CONCLUSIONS 



We have built a complete computer model of the Oklo reactor core RZ2 . With the aid of present-day computational 
codes we have calculated in all detail the core parameters. The simulations were done for three fresh cores of different 
contents of uranium and water. We have also calculated the neutron flux and its spatial and energy distributions. 
For the three cores we have estimated the temperature and void effects in the reactor. As expected, the neutron 
reactor spectrum is significantly different from the ideal Maxwell distribution that had been used by other authors 
to determine the cross section of ^If^™- The reactor cross section and the curves of its dependence on the shift of 
the resonance position AEr (as a result of a possible change of fundamental constants) differ appreciably from earlier 
results. The effect of an influence of the reactor spectrum on the cross section of ^52^™ be considered to be 
firmly established. We have studied the limits of the variation of this effect depending on the initial composition and 
the size of the core. The fresh bare core RZ2 is critical for Tc = (725 ± 55)K. At these temperatures the curves of 
'^■y,Sm{Tc, /^Er) lie appreciably lower than for a Maxwell distribution. Possible values of AEr lie in the range of —73 
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TABLE XIL Limits on the rate of a variation, based on the data of ^g^Sm contents in the Oklo reactor 
{z = 1.8 • 10^year/13.7(2) • 10^year= 0.131(2)); to = 13.7(2) [l^l and obtained by other methods. 





Lab 


Authors, year 


A_Br, 5a/a, (5a/a 


Comments 
core, spectrum, Tc 


1 


LNPI, Gatchina 
Russia 


A.Shlyakhter, 
1976 [8] 


AEr < 20 meV 
la/ a < 0.5 ■ 10"" year"^ 


RZ2 , MaxweU, 300 K 


2 


LNPI, 

1 -Ici^"/^ Vi in CI 

Russia 


Yu. Petrov, 
1 Q77 \^ 01 


AEr < 50 meV 


RZ2 , Maxwell, 300K 


3 


Princeton, 
USA 


T.D amour and 

F.Dyson, 
1 QQfi fn 1 21 


Ai;^ < 90 meV 
5Q/a < 5.0 • 10"^^ year"^ 


RZ2 , Maxwell, (450-1000)K 


4 


Univ.Tokio, 

Tokyo, 
J&pcLn 


Ya.Fujii 
et al.. 


AEr < 20 meV 
la/a < 1.0 • 10"" year"^ 


RZIO, MaxweU, (470-670) K 


5 


LANL, 
USA 


S.Lamoreaux and 

J . ± (Ji g"! bUli 5 

2004 [g^ 


AEr < -45_+^ meV 

ua / a ^ O.O -LU yeai 


RZIO, Maxwell+Fermi, 

UUU Iv 




This paper. 

Fig. uni 


AB,. < 4 meV 
5Q/a < 0.2 ■ 10~^^ year"^ 


RZIO, Reactor spectrum 

for LT core, 
600 K 


6 


PNPI, 

Gatchina, 
Russia 


This paper 


AEr < 62 meV 
la/a < 3.1 ■ 10~" year"^ 


RZ2 , Reactor spectrum 

for fresh core, 
725 ± 55 K 


Cosmophysical and laboratory data 


7 


lUCAA, Pune, 
India 


H.Chand, 

R.Srianand et al., 
2004 


Sa/a < (-60 ±60) • 10^** 
S^/a < 12 • 10"^'' year-i 


Cosmophysical 
multidoublet method 


8 


Observ. 
de Paris, 
France 


S.Bize et al. 
2004 


5a/a < (-5 ± 53) • 10"^'^ year~^ 


Method of 
Atomic Fountains 



meV < AEr < 62 meV. These limits are 1.5 times more accurate than those of Dyson and Damour. For the rate of 

change of the e.m. constant we find —3.7 • 10^^^ year^^ < 5a/a < +3.1 • 10^^^ year^^. Within these limits we have 

Sa/a = 0. The analysis of all previous studies shows that none of them has reliably shown up a deviation from zero 
of the rate of change of the e.m. constant a. Because of difficulties with the detailed calculation of the burning up in 
large reactors, which require accumulation of huge statistics, we have not determined the effect of the burn-up on the 
neutron spectrum and on the ^If^m cross section. Calculations of the influence of burn-up on the temperature of the 
active core and on the neutron spectrum are in progress. 
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